Method for Determining Azimuth and Elevation Angles of Arrival of Coherent Sources

ABSTRACT

A method for jointly determining the azimuth angle θ and the elevation angle Δ of the wave vectors of P waves in a system comprising an array of sensors, a number of waves out of the P waves being propagated along coherent or substantially coherent paths between a source and said sensors, includes at least the following steps: selecting a subset of sensors from said sensors to form a linear subarray of sensors; applying, to the signals from the chosen subarray, an algorithm according to a single dimension to decorrelate the sources of the P waves; determining a first component w of said wave vectors by applying, to the signals observed on the sensors of the chosen subarray, a goniometry algorithm according to the single dimension w; determining a second component u of said wave vectors by applying a goniometry algorithm according to the single dimension u to the signals from the entire array of sensors; determining θ and Δ from w and u.

The present invention relates to a method for determining azimuth andelevation angles of arrival of coherent sources. It is used, forexample, in all locating systems in an urban context in which thepropagation channel is disturbed by a large number of obstacles such asbuildings. Generally, it can be used to locate transmitters, for examplea cell phone, in a difficult propagation context, whether it be an urbanenvironment, a semi-urban environment—for example an airport site—, theinterior of a building or under snow, in the case of an avalanche. Theinvention can also be used in medical imaging methods, to locate tumorsor sources of epilepsy, the tissues of the human body being the originof multiple wave paths. It also applies to sounding systems for miningand oil surveys in the seismic domain, in which the aim is to estimateangles of arrival with multiple paths in the complex propagation mediumof the earth's crust.

The invention is situated in the technical field of antenna processingwhich processes the signals from a number of transmitting sources basedon a multiple-sensor reception system. More specifically, the inventionrelates to the field of goniometry which consists in estimating theangles of arrival of the sources.

In an electromagnetic context, the sensors are antennas and theradiofrequency sources are propagated according to a polarization thatis dependent on the transmitting antenna. In an acoustic context, thesensors are microphones and the sources are sounds.

FIG. 1 shows that an antenna processing system comprises an array 102 ofsensors receiving sources with different angles of arrival θ_(mp). Theindividual sensors 101 of the array receive the sources E1, E2 with aphase and an amplitude that are dependent in particular on their anglesof incidence and on the position of the sensors. The angles of incidenceare parameterized in one dimension (1D) by the azimuth θ_(m) and in twodimensions (2D) by the azimuth θ_(m) and the elevation Δ_(m).

According to FIG. 2, a 1D goniometry is defined by techniques whichestimate only the azimuth by assuming that the waves of the sources arepropagated in the plane 201 of the array of sensors. When the goniometrytechnique jointly estimates the azimuth and elevation of a source, it isan issue of 2D goniometry.

The aim of the antenna processing techniques is notably to exploit thespatial diversity generated by the multiple-antenna reception of theincidence signals, in other words, to use the position of the antennasof the array to better use the differences in incidence and distance ofthe sources.

FIG. 3 illustrates an application to goniometry in the presence ofmultiple paths. The m-th source 301 is propagated along P paths 311,312, 313 of incidences θ_(mp) (1≦p≦P) which are provoked by P−1obstacles 320 in the radiofrequency environment. The problem dealt within the method according to the invention is, notably, to perform a 2Dgoniometry for coherent paths in which the propagation time deviationbetween the direct path and a secondary path is very low.

One known method for doing the goniometry is the MUSIC algorithm [1].However, this algorithm does not make it possible to estimate theincidences of the sources in the presence of coherent paths.

The algorithms that make it possible to process the case of coherentsources are the maximum likelihood algorithms [2][3] which areapplicable to arrays of sensors with any geometry. However, thesetechniques require the calculation of a multidimensional criteria ofwhich the number of dimensions depends on the number of paths and on thenumber of incidence parameters for each path. More particularly, in thepresence of K paths, the criterion has 2K dimensions for a 2D goniometryin order to jointly estimate all the incidences (Θ₁, . . . , Θ_(K)). Itshould be noted that, even in the presence of a number K′ of coherentpaths that is less than the total number K of paths, the calculation ofthe maximum likelihood criterion still has 2K dimensions. In order toreduce the number of dimensions of the criterion to 2K′, one alternativeis to apply the coherent MUSIC method [4]. However, the coherent MUSICalgorithm [4] requires a high number of sensors and very significantcomputation resources.

Another alternative for reducing the computation cost is to implementspatial smoothing or forward-backward techniques [5][6], thesetechniques requiring particular array geometries. In practice, spatialsmoothing is applicable when the array is broken down into subarrayshaving the same geometry (e.g.: evenly-spaced linear array or array on aregular 2D grid). The forward-backward algorithm requires an array witha center of symmetry. These techniques are highly restrictive in termsof geometry of the array of sensors, especially for a 2D goniometry, inwhich the constraint of symmetry or of translated identical subarrays isdifficult to satisfy.

For 1D goniometry, spatial smoothing techniques have been considered onany arrays [6][7]. For this, the array of sensors is interpolatedaccording to the goniometry adapted to spatial smoothing (orforward-backward). In [6] the interpolation technique addresses only asingle angular segment and in [7] the algorithm is adapted to the casesof a number of angular segments for the interpolation. However, thiskind of technique is difficult to adapt to the case of 2D goniometry.

One aim of the invention is to propose a method for determining, from anarray of sensors, the direction of arrival in azimuth and elevation ofcoherent signals with a reduced computation cost and by limiting as faras possible the geometrical constraints to be imposed on the array ofsensors. To this end, the subject of the invention is a method forjointly determining the azimuth angle θ and the elevation angle Δ of thewave vectors of P waves in a system comprising an array of sensors, anumber of waves out of the P waves being propagated along coherent orsubstantially coherent paths between a source and said sensors, themethod being characterized in that it comprises at least the followingsteps:

selecting a subset of sensors from said sensors to form a linearsubarray of sensors;

applying, to the signals from the chosen subarray, an algorithmaccording to a single dimension to decorrelate the sources of the Pwaves;

determining a first component w of said wave vectors by applying, to thesignals observed on the sensors of the chosen subarray, a goniometryalgorithm according to the single dimension w;

determining a second component u of said wave vectors by applying agoniometry algorithm according to the single dimension u to the signalsfrom the entire array of sensors;

determining θ and Δ from w and u.

The method according to the invention makes it possible to reduce thecomplexity of a two-dimensional goniometry problem by dividing it intotwo phases with a single dimension: a first phase for estimating aprojection value w of the wave vector, then a phase for estimating avalue of the other component u of the wave vector.

According to one implementation of the method according to theinvention,

the sensors forming the subarray are chosen such that at least a portionof the subarray is unchanging by translation;

a spatial smoothing algorithm is applied to decorrelate the sources ofthe P waves.

The expression “unchanging by translation” should be understood to meanthat there is at least one subset of the subarray whose sensors arearranged as if they were the result of a translation of another subsetof sensors of said subarray.

According to one implementation of the method according to theinvention,

the sensors forming the subarray are chosen such that at least a portionof the subarray includes a center of symmetry;

a forward-backward algorithm is applied to decorrelate the sources ofthe P waves.

The expression “center of symmetry of the subarray” should be understoodto mean that, for each sensor, there is a sensor placed symmetricallyrelative to said center.

According to one implementation of the method according to theinvention, the determined first component w of the wave vectors is theprojection, on the axis formed by the linear subarray, of the projectionof the wave vectors on the plane formed by the array of sensors. Inother words, for each path p, w_(p)=cos(θ_(p)−α).cos(Δ_(p)), α being theazimuth angle according to which the axis formed by the linear subarrayis oriented.

According to one implementation of the method according to theinvention, the method includes at least the following steps:

calculating the covariance matrix R_(x) on the entire array of sensors;

extracting from R_(x), the covariance matrix R_(x)′ corresponding to thechosen linear subarray;

applying a source decorrelation algorithm to R_(x)′;

estimating, for each path p, the values of the first componentw_(p)=cos(θ_(p)−α).cos(Δ_(p)) by applying a 1D goniometry algorithm tothe decorrelated matrix R_(x)′, α being the azimuth orientation angle ofthe axis formed by the linear subarray;

estimating the values of the second componentu_(p)=cos(θ_(p)).cos(Δ_(p)), for each path p, by applying a 1Dgoniometry algorithm to the matrix R_(x);

determining, from the values of the pairs (w_(p), u_(p)), the values ofthe azimuth-elevation pairs (θ_(p), Δ_(p)).

According to one implementation of the method according to theinvention, the goniometry algorithm used to determine the firstcomponent w of each wave vector P is the MUSIC algorithm, the criterionJ_(MUSIC) to be minimized to determine said component w being equal to

${{J_{MUSIC}(w)} = \frac{{a^{l}(w)}^{H}\Pi_{b}^{l}{a^{l}(w)}}{{a^{l}(w)}^{H}{a^{l}(w)}}},$

in which Π_(b) ^(l) is the noise projector extracted from thedecorrelated covariance matrix R_(x)′ corresponding to the chosen linearsubarray, and a(w)^(l) represents the response of the chosen subarray tothe incident waves P.

According to one implementation of the method according to theinvention, the goniometry algorithm used to determine the secondcomponent u is the coherent MUSIC algorithm in a single dimension, thecriterion to be minimized being:

${J_{{coherent}\mspace{14mu} {MUSIC}}\left( \underset{\_}{\Theta} \right)} = \frac{\det \left( {{D\left( \underset{\_}{\Theta} \right)}^{H}\Pi_{b}{D\left( \underset{\_}{\Theta} \right)}} \right)}{\det \left( {{D\left( \underset{\_}{\Theta} \right)}^{H}{D\left( \underset{\_}{\Theta} \right)}} \right)}$

in which Θ={f₁(u₁) . . . f_(K) _(max) (u_(K) _(max) )}, with

${{f_{p}(u)} = {f\left( {u,\frac{w_{p} - {u\; {\cos (\alpha)}}}{\sin \; (\alpha)}} \right)}},$

α being the azimuth orientation angle of the axis formed by the linearsubarray, D(Θ) being a vector equal to [a(Θ ₁) . . . a(Θ _(Kmax))], a(Θ_(i)) being the response of the array of sensors to the path of index i,K_(max) being the maximum number of coherent paths.

According to one implementation of the method according to theinvention, the goniometry algorithm used to determine the secondcomponent is the maximum likelihood algorithm.

According to one implementation of the method according to theinvention, the array is disturbed by mutual coupling of known matrix Z,and the method includes a step for eliminating the coupling executedprior to the steps for estimating the values of the components w and u,said step for eliminating the coupling determining a covariance matrixthat is cleaned of noise by applying the following processing to thecovariance matrix: Z⁻¹(R_(x)−σ²|)Z^(−1H)σ² being the estimated noiselevel.

According to one implementation of the method according to theinvention, the determination of the pairs of values (θ_(p), Δ_(p)) fromthe values of pairs (w_(p), u_(p)) is performed as follows:

$\quad\left\{ \begin{matrix}{\theta_{p} = {{angle}\left( {u_{p} + {j\; v_{p}}} \right)}} \\{\Delta_{p} = {\cos^{- 1}\left( {{u_{p} + {j\; v_{p}}}} \right)}}\end{matrix} \right.$

in which

$v_{p} = {\frac{w_{p}u_{p}\cos \; (\alpha)}{\sin \; (\alpha)}.}$

Other characteristics will emerge from reading the following detaileddescription given as a nonlimiting example and in light of the appendeddrawings which represent:

FIG. 1, an example of signals transmitted by a transmitter and beingpropagated to an array of sensors,

FIG. 2, the representation of the incidence of a source on a plane ofsensors,

FIG. 3, an illustration of the propagation of signals by multiple paths,

FIG. 4, an example of arrays of sensors of position (x_(n),y_(n)),

FIG. 5, an example of an array of sensors consisting of two subarraysthat are unchanging by translation,

FIG. 6, an example of an array of sensors consisting of two subarraysmaking it possible to decorrelate two paths that are coherent in azimuthand elevation,

FIG. 7, a linear array of sensors in which the spatial smoothing makesit possible to decorrelate two paths that are coherent in azimuth,

FIG. 8, an array of sensors having a center of symmetry at O,

FIG. 9, two linear arrays with five sensors, making it possible todecorrelate two paths that are coherent in azimuth for, respectively,spatial smoothing and forward-backward,

FIG. 10, a first example of an array of sensors containing a linearsubarray of sensors and compatible with the method according to theinvention,

FIG. 11, a second example of an array of sensors containing a linearsubarray of sensors and compatible with the method according to theinvention.

Before detailing an exemplary implementation of the method according tothe invention, some reminders concerning the modeling of the outputsignal from an array of sensors are given.

With M sources in which the m-th source contains P_(m) multiple paths,the output signal from the array of sensors is written as follows:

$\begin{matrix}{{x(t)} = {\begin{bmatrix}{x_{1}(t)} \\\vdots \\{x_{N}(t)}\end{bmatrix} = {{\sum\limits_{m = 1}^{M}\; {\sum\limits_{p = 1}^{P_{m}}\; {\rho_{mp}{a\left( {\underset{\_}{\Theta}}_{mp} \right)}{s_{m}\left( {t - \tau_{mp}} \right)}}}} + {{n(t)}.}}}} & (1)\end{matrix}$

in which x_(n)(t) is the signal at the output of the n-th sensor, N thenumber of sensors, n(t) is the additive noise, a(Θ) is the response ofthe array of sensors to a source of direction Θ=(θ, Δ), θ is theazimuth, Δ the elevation and ρ_(mp), θ_(mp), τ_(mp) are respectively theattenuation, the direction and the delay of the p-th paths of the m-thsource. The vector a(Θ) which is also called directing vector depends onthe positions (x_(n), y_(n)) of the sensors 401, 402, 403, 404, 405 (seeFIG. 4) and is written

$\begin{matrix}{{a\left( \underset{\_}{\Theta} \right)} = {{Z\begin{bmatrix}{a_{1}\left( \underset{\_}{\Theta} \right)} \\\vdots \\{a_{N}\left( \underset{\_}{\Theta} \right)}\end{bmatrix}}\mspace{14mu} {with}\mspace{14mu} \left\{ \begin{matrix}{{a_{n}\left( \underset{\_}{\Theta} \right)} = {\exp \left( {j\frac{2\pi}{\lambda}\left( {{x_{n}u} + {y_{n}v}} \right)} \right)}} \\{u = {{\cos (\theta)}{\cos (\Delta)}}} \\{v = {{\sin (\theta)}{{\cos (\Delta)}.}}}\end{matrix} \right.}} & (2)\end{matrix}$

in which Z is the coupling matrix, λ is the wavelength and (u,v) are thecoordinates of the wave vector in the plane of the antenna.

With coherent paths in which the delays of the paths satisfy τ_(m1)= . .. =τ_(mPm,) the signal model of the equation (1) becomes

$\begin{matrix}{{{x(t)} = {{\sum\limits_{m = 1}^{M}\; {{a\left( {{\underset{\_}{\Theta}}_{m},\rho_{m},P_{m}} \right)}{s_{m}(t)}}} + {n(t)}}}{with}{{a\left( {{\underset{\_}{\Theta}}_{m},\rho_{m},P_{m}} \right)} = {\sum\limits_{p = 1}^{P_{m}}\; {\rho_{mp}{{a\left( {\underset{\_}{\Theta}}_{mp} \right)}.}}}}} & (3)\end{matrix}$

in which a(Θ _(m), ρ_(m), P_(m)) is the response of the array of sensorsto the m-th source, Θ _(m)=[Θ _(m1) . . . Θ _(mP) _(m) ]^(T) andρ_(m)=[ρ_(m1) . . . ρ_(mP) _(m) ]^(T). The directing vector of thesource is no longer a(Θ _(m1)) but a composite directing vector a(Θ_(m),ρ_(m),P_(m)) dependant on a greater number of parameters.

More generally, with K groups of coherent paths, the signal is written:

$\begin{matrix}{{{x(t)} = {{\sum\limits_{k = 1}^{K}\; {{a\left( {{\underset{\_}{\Theta}}_{k},\rho_{k},K_{k}} \right)}{s_{k}(t)}}} + {n(t)}}}{{in}\mspace{14mu} {which}}{K_{\max} = {\max\limits_{k}{\left\{ K_{k} \right\}.}}}} & (4)\end{matrix}$

To enable the reader to better understand the method according to theinvention, the processing of the coherent sources in azimuth andelevation in the state of the art is explained hereinbelow.

A first coherent MUSIC algorithm [4] is first described. The MUSICalgorithm [1] is a high-resolution method based on the breakdown intospecific elements of the covariance matrix R_(x)=E[x(t) x(t)^(H)] of themultiple-sensor signal x(t), in which E[.] is the mathematicalexpectation. The expression of the matrix R_(X) is as follows accordingto (4):

$\begin{matrix}{{R_{x} = {{{AR}_{s}A^{H}} + {\sigma^{2}I_{N}}}}{{with}\mspace{14mu} \left\{ \begin{matrix}{R_{s} = {E\left\lbrack {{s(t)}{s(t)}^{H}} \right\rbrack}} \\{{E\left\lbrack {{n(t)}{n(t)}^{H}} \right\rbrack} = {\sigma^{2}I_{N}}} \\{A = \left\lbrack {{a\left( {{\underset{\_}{\Theta}}_{1},\rho_{1},K_{1}} \right)}\mspace{14mu} \ldots \mspace{14mu} {a\left( {{\underset{\_}{\Theta}}_{k},\rho_{k},K_{k}} \right)}} \right\rbrack} \\{{s(t)} = {\left\lbrack {{s_{1}(t)}\mspace{14mu} \ldots \mspace{14mu} {s_{K}(t)}} \right\rbrack^{T}.}}\end{matrix} \right.}} & (5)\end{matrix}$

With K groups of coherent paths, the rank of the matrix R_(x) is K. Inthese conditions, the K specific vectors e_(k) (1≦k≦K) associated withthe K highest specific values λ_(k) of R_(x) satisfy

$\begin{matrix}{e_{k} = {\sum\limits_{i = 1}^{K}{\alpha_{ik}{a\left( {{\underset{\_}{\Theta}}_{i},\rho_{i},K_{i}} \right)}\mspace{14mu} {for}\mspace{14mu} {\left( {1 \leq k \leq K} \right).}}}} & (6)\end{matrix}$

The N-K specific vectors e_(i)(K+1≦l≦N) associated with the lowestspecific values of R_(x) are orthogonal to the vectors e_(k) (1≦k≦K)ofthe expression (6) and define the noise space. Since the vectors e_(i)and e_(k) are orthogonal, the directing vectors a(Θ _(i), ρ_(i), K_(i))are orthogonal to the noise vectors e_(i). In these conditions, the Kminima (Θ _(k), ρ_(k),K_(k)) of the following MUSIC criterion

$\begin{matrix}{{{J_{MUSIC}\left( {\underset{\_}{\Theta},\rho} \right)} = \frac{{a\left( {\underset{\_}{\Theta},\rho,K_{\max}} \right)}^{H}\Pi_{b}{a\left( {\underset{\_}{\Theta},\rho,K_{\max}} \right)}}{{a\left( {\underset{\_}{\Theta},\rho,K_{\max}} \right)}^{H}{a\left( {\underset{\_}{\Theta},\rho,K_{\max}} \right)}}}{{{with}\mspace{14mu} \Pi_{b}} = {\sum\limits_{i = {K + 1}}^{N}\; {e_{i}{e_{i}^{H}.}}}}} & (7)\end{matrix}$

make it possible to give the directions Θ _(k) of each path. However,the cost of calculating the criterion of the equation (7) is very highbecause it depends on the incidence of K_(max) coherent paths and theirrelative amplitudes: (Θ, ρ).

The coherent MUSIC method described in [4] is designed to reduce thenumber of parameters for searching for the MUSIC criterion. For this,the vector a(Θ _(m), ρ_(m), P_(m)) of equation (3) is rewritten asfollows:

$\begin{matrix}{{{a\left( {\underset{\_}{\Theta},\rho,K_{\max}} \right)} = {{D\left( \underset{\_}{\Theta} \right)}\rho}}{{with}\mspace{14mu} \left\{ \begin{matrix}{{D\left( \underset{\_}{\Theta} \right)} = \left\lbrack {{a\left( {\underset{\_}{\Theta}}_{1} \right)}\mspace{14mu} \ldots \mspace{14mu} {a\left( {\underset{\_}{\Theta}}_{K_{\max}} \right)}} \right\rbrack} \\{\rho = \left\lbrack {\rho_{1}\mspace{14mu} \ldots \mspace{14mu} \rho_{K_{\max}}} \right\rbrack^{T}} \\{\underset{\_}{\Theta} = {\left\{ {{\underset{\_}{\Theta}}_{1}\mspace{14mu} \ldots \mspace{14mu} {\underset{\_}{\Theta}}_{K_{\max}}} \right\}.}}\end{matrix} \right.}} & (8)\end{matrix}$

In these conditions, the criterion of equation (7) is reduced to thefollowing expression:

$\begin{matrix}{{J_{{coherent}\mspace{14mu} {MUSIC}}\left( \underset{\_}{\Theta} \right)} = {\frac{\det \left( {{D\left( \underset{\_}{\Theta} \right)}^{H}\Pi_{b}{D\left( \underset{\_}{\Theta} \right)}} \right)}{\det \left( {{D\left( \underset{\_}{\Theta} \right)}^{H}{D\left( \underset{\_}{\Theta} \right)}} \right)}.}} & (9)\end{matrix}$

in which det(M) is the determinant of the matrix M. The number ofdimensions of the criterion is then reduced to 2K_(max) parameters of Θin which K_(max) is the maximum number of coherent paths. Consequently,the K minima of the criterion J_(MUSIC-Coherent) (Θ) gives thedirections Θ _(k)={Θ _(k1) . . . Θ _(kK) _(max) } of the paths of eachgroup of coherent paths for 1≦k≦K . Thus, with K_(max)=2 coherent paths,the coherent MUSIC criterion still has four dimensions for anazimuth-elevation goniometry. More generally, the coherent MUSIC methodentails calculating a criterion (9) having 2 K_(max) dimensions.However, the method makes no assumption as to the geometry of the arraybecause it entails no constraint on the expression of the directingvector a(Θ).

Alternative methods that are also known for processing coherent sourcesare the spatial smoothing [5][6] and forward-backward [5] techniques.These methods make it possible to decorrelate the sources by performinga simple preprocessing on the covariance matrix of the received signals.It is then possible to apply a goniometry algorithm such as MUSIC to thenew covariance matrix. These techniques derive from the field ofspectral analysis whose objective is to model the frequency spectrum ofa signal.

The spatial smoothing techniques [5][6] are applicable to an array ofsensors consisting of subarrays 501, 502 that are unchanging bytranslation as illustrated in FIG. 5. With P paths (coherent or not),the expression (1) of the observation vector can be rewritten

$\begin{matrix}{{x(t)} = {{{\sum\limits_{p = 1}^{P}\; {\rho_{p}{a\left( {\underset{\_}{\Theta}}_{p} \right)}{s_{p}(t)}}} + {n(t)}} = {{{As}(t)} + {n(t)}}}} & (10)\end{matrix}$

with A=[a(Θ ₁) . . . a(Θ _(P))].

The expression of the signal received on the i-th subarray is thenwritten:

$\begin{matrix}{{x^{i}(t)} = {{P^{i}{x(t)}} = {{{\sum\limits_{p = 1}^{P}\; {\rho_{p}{a^{i}\left( {\underset{\_}{\Theta}}_{p} \right)}{s_{p}(t)}}} + {n(t)}} = {{A^{i}\mspace{14mu} {s(t)}} + {n(t)}}}}} & (11)\end{matrix}$

in which A^(i)=[a^(i)(Θ ₁) . . . a^(i)(Θ _(P))], P^(i) being a matrixconsisting of 0 and 1 making it possible to select the signal of thei-th subarray for which the directing vector a^(i)(Θ) satisfies thefollowing relationship:

a ^(i)(Θ)=P ^(i) a(Θ)=α^(i)(Θ)a¹(Θ)  (12)

Remembering that the incidence Θ=(θ,Δ) depends on the two parameters θand Δ.

According to (11)(12), the mixing matrix A^(i) of the i-th subarraysatisfies

A ^(i) =P ^(i) A=A ¹Φ_(i) with Φ_(i)=diag{α^(i)(Θ ₁) . . . α^(i)(Θ_(P))}  (13)

According to (11)(13) the covariance matrix R_(x)^(i)=E[x(t)^(i)x(t)^(iH)] has the following expression:

R _(x) ^(i) =A ¹Φ_(i) R _(s)Φ_(i) *A ^(1H)+σ² I _(N) in which R _(s)=E[s(t) s(t)^(H)]  (14)

Consequently, an alternative to the spatial smoothing techniquesconsists in applying a MUSIC-type algorithm to the following covariancematrix:

$\begin{matrix}{{R_{x}^{SM}{\sum\limits_{i = 1}^{I}\; R_{x}^{i}}} = {\sum\limits_{i = 1}^{I}\; {P^{i}{R_{x}\left( P^{i} \right)}^{H}}}} & (15)\end{matrix}$

in which R_(x)=E[x(t) x(t)^(H)]. The aim of this procedure is to obtaina matrix R_(x) ^(SM) that has a rank higher than the R_(x) ^(i) withoutdestroying the structure of the signal space generated by A¹. Inpractice, this technique makes it possible to decorrelate a maximumof/coherent paths because

$\begin{matrix}{{R_{x}^{SM} = {{A^{1}R_{s}^{SM}A^{1\; H}} + {\sigma^{2}I_{N^{\prime}}\mspace{14mu} {in}\mspace{14mu} {which}}}}{R_{s}^{SM} = {\sum\limits_{i = 1}^{I}{\Phi_{i}R_{s}\Phi_{i}^{*}}}}} & (16)\end{matrix}$

and thus, rank{R_(s)}≦rank{R_(s) ^(SM)}≦min(rank{R_(s)}/, Σ_(m=1)^(M)P_(m)).

The forward-backward [5] smoothing technique requires an array ofsensors that has a center of symmetry at O as indicated in FIG. 8. Inthese conditions, the directing vector has the following structure

$\begin{matrix}{{a\left( \underset{\_}{\Theta} \right)} = {{\beta \left( \underset{\_}{\Theta} \right)}\begin{bmatrix}{b\left( \underset{\_}{\Theta} \right)} \\{b\left( \underset{\_}{\Theta} \right)}^{*}\end{bmatrix}}} & (17)\end{matrix}$

in which, according to FIG. 8, b(Θ) is the directing vector of thesubarray of coordinates (x_(n)-x₀,_(yn)-y₀) and b(Θ)* is the directingvector of the subarray of coordinates (-x_(n)-x₀,-y_(n)-y₀), bearing inmind that (x₀, y₀) are the coordinates of the center of symmetry O.Consequently, the directing vector of the expression (17) satisfies thefollowing relationship:

Πa(Θ)*=β(Θ)* a(Θ)  (18)

in which Π is a permutation matrix consisting of 0 and 1. Theforward-backward smoothing technique consists in applying a goniometryalgorithm such as MUSIC to the following covariance matrix

R _(x) ^(FB) =R _(x) +ΠR _(x)*Π^(T)  (19)

bearing in mind that

R _(x) ^(FB) =AR _(s) ^(FB) A ^(H)+σ² I _(N) in which R _(s) ^(FB) =R_(s)+Φ_(FB) R _(s)Φ_(FB)*  (20)

The technique makes it possible to decorrelate two coherent pathsbecause rank{R_(s)}≦rank{R_(s) ^(SM)}≦min(2rank{R_(s)}, Σ_(m=1)^(M)P_(m)) with

Φ_(FB)=diag{β(Θ ₁) . . . β(Θ _(P))}  (21)

The spatial smoothing and forward-backward techniques can be combined toincrease the capacity for decorrelation into number of paths. Thesesmoothing techniques make it possible to process the coherent sourceswith a computation power that is very similar to the application of asingle goniometry algorithm such as MUSIC.

When the array of sensors is disturbed by mutual coupling in which thedirecting vector is written

a(Θ)=Z

(Θ)  (22)

and in which the directing vector

(Θ) satisfies one of the properties of the equations (12)(18), thespatial smoothing techniques are applicable [7]. The mixing matrix A ofthe equation (10) is then written

A=Z

with

=[

(Θ ₁) . . .

(Θ _(P))]  (23)

Consequently, the covariance matrix R_(x)=E[x(t) x(t)^(H)] is written asfollows:

R _(x) =Z(

R _(s)

^(H))Z ^(H)+σ² I _(N)  (24)

Bearing in mind that

^(i)=P^(i)

=

¹Φ_(i) (where that Π

*=

Φ_(FB)), the following steps make it possible to apply a spatialsmoothing or forward-backward technique in the presence of mutualcoupling:

Step No. L.1: Break down the covariance matrix R_(x)=E[x(t) x(t)^(H)]into specific elements such that:

R _(x) =E _(s)Λ_(s) E _(s) ^(H) +E _(b)Λ_(b) E _(b) ^(H)  (25)

in which E_(s) and E_(b) are the matrices of the specific vectorsrespectively associated with the signal space and the noise spaceaccording to MUSIC0 and in which Λ_(s) and Λ_(b) are diagonal matricesrespectively consisting of the specific values of the signal space andof the specific values of the noise space.

Step No. L.2 : Extract the non-noise-affected covariance matrix Z(

R_(s)

^(H))Z^(H) by performing:

$R_{y} = {{R_{x} - {\frac{{trace}\left( \Lambda_{b} \right)}{N - K}I_{N}}} = {{Z\left( {\overset{\Cap}{A}R_{s}{\overset{\Cap}{A}}^{H}} \right)}Z^{H}}}$

in which K is the dimension of the signal space such that K≦P.

Step No. L.3 (Spatial smoothing): Apply the MUSIC algorithm to thefollowing covariance matrix R_(x) ^(SM):

$R_{x}^{SM} = {\sum\limits_{i = 1}^{I}\; {{P^{i}\left( {Z^{- 1}{R_{y}\left( Z^{- 1} \right)}^{H}} \right)}\left( P^{i} \right)^{H}}}$

Step No. L.3 (Forward-Backward): Apply the MUSIC algorithm to thefollowing covariance matrix R_(x) ^(FB):

R _(x) ^(FB)=(Z ⁻¹ R _(y)(Z ⁻¹)^(H))+Π(Z ⁻¹ R _(y)(Z ⁻¹)^(H))*Π^(T)

If the directing vector array

(Θ) permits, the two smoothing techniques of steps No. L.3 can becombined.

The spatial smoothing techniques are applicable with mutual coupling.However, this imposes very strong constraints on the geometry of theindividual array which have the drawback of requiring a very largenumber of sensors. In the following example, we will evaluate theminimum number of sensors to process the case of two sources coherent inazimuth-elevation. For this, it is necessary for:

Constraint C1: The number of sensors of each subarray to be at leastequal to N^(i)=4. In practice, because of ambiguities, an array of Nsensors makes it possible at most to estimate the direction of arrivalof N/2 sources.

Constraint C2: The number of subarrays to be at least equal to 2.

Constraint C3: The subarrays to be planar (not linear) in order to beable to perform an azimuth-elevation goniometry.

FIG. 6 shows that an array consisting of two subarrays 601, 602 of foursensors contains at least seven sensors. This array also has thedrawback of being weakly open (or has little spatial bulk) because thesubarrays with four sensors must be unambiguous. Since the subarraysconsist of four sensors, this ambiguity constraint requires a spacingbetween sensors less than λ/2. In practice, the more open a array is,the more accurate the estimation of the angles of arrival is with abetter robustness to calibration errors. For the case where the desireis to perform an azimuth goniometry only, the constraint C3 no longerapplies and the array making it possible to perform a goniometry on twocoherent sources consisting of two subarrays of four sensors is anevenly-spaced linear array with five sensors. Each subarray is then anevenly-spaced linear array with four sensors. FIG. 7 shows that thelinear subarray 701 which allows for the azimuth goniometry of twocoherent paths has the following differences compared to the array 702which makes it possible to do so in azimuth and elevation: on the onehand, it consists of fewer sensors: five instead of seven, and on theother hand, it has a greater bulk: 4d instead of 3d, bearing in mindthat d is a distance less than λ/2, λ being the wavelength of thetransmitted signals.

For the forward-backward technique requiring an array with a center ofsymmetry as illustrated in FIG. 8, it is possible to note that, for thespatial smoothing: the decorrelation of two coherent paths for anazimuth-elevation goniometry requires an array of sensors having moresensors and less aperture than the array making it possible to performan azimuth goniometry only. For an azimuth goniometry, a linear array,not necessarily evenly-spaced, is sufficient.

FIG. 9 shows that the forward-backward technique makes it possible,compared to the spatial smoothing technique, to perform an unambiguousgoniometry on two coherent paths with an array 901 having a greateraperture (10d instead of 4d for an array 902 used for spatialsmoothing). The forward-backward technique has the advantage of notimposing any geometry constraint on half the array. The other half ofthe array is symmetrical to the 1st half.

The method according to the invention described combines the coherentMUSIC method with a forward-backward technique and/or a spatialsmoothing technique. Given the advantages and drawbacks of the smoothingtechniques and of the coherent MUSIC algorithm described above, themethod envisages using an array of sensors containing a linear subarrayon which a spatial smoothing and/or forward-backward technique can beenvisaged. FIG. 10 shows such an array 1001 with a linear subarray 1002having an orientation a relative to the x axis. More specifically, themethod according to the invention can use the array 1101 of FIG. 11 inwhich the angle α=90° and the linear subarray 1102 consists of 3evenly-spaced sensors on which a forward-backward technique can beapplied.

The coordinates (x_(n) ^(l), y_(n) ^(l)) of the n-th sensor of thelinear subarray then have the following expression:

$\begin{matrix}\left\{ {{\begin{matrix}{x_{n}^{l} = {\rho_{n}{\cos (\alpha)}}} \\{y_{n}^{l} = {\rho_{n}{\sin (\alpha)}}}\end{matrix}{for}\mspace{14mu} 1} \leq n \leq N^{l}} \right. & (26)\end{matrix}$

in which N^(l) is the number of sensors of the linear subarray. In theabsence of coupling and according to (2), the directing vector a^(l)(Θ)associated with the linear subarray is written

$\begin{matrix}{{a^{l}\left( \underset{\_}{\Theta} \right)} = {{\begin{bmatrix}{a_{1}^{l}\left( \underset{\_}{\Theta} \right)} \\\vdots \\{a_{N^{l}}^{l}\left( \underset{\_}{\Theta} \right)}\end{bmatrix}\mspace{14mu} {with}\mspace{14mu} {a_{n}^{l}\left( \underset{\_}{\Theta} \right)}} = {\exp \left( {{j2\pi}\frac{\rho_{n}}{\lambda}{\cos \left( {\theta - \alpha} \right)}{\cos (\Delta)}} \right)}}} & (27)\end{matrix}$

The vector a^(l)(Θ) then depends on a single parameter w=cos(θ−α)cos(Δ)as follows:

$\begin{matrix}{{a^{l}\left( \underset{\_}{\Theta} \right)} = {{a^{l}(w)} = {{\begin{bmatrix}z^{\rho_{1}} \\\vdots \\z^{\rho_{N^{l}}}\end{bmatrix}\mspace{14mu} {with}\mspace{14mu} z} = {\exp \left( {j\frac{2\pi}{\lambda}w} \right)}}}} & (28)\end{matrix}$

x^(l)(t) is used to denote the signal at the output of the linearsubarray and P_(roj) the matrix consisting of 0 and 1 that can be usedto extract the signals from the linear subarray such that

x ^(l)(t)=P _(roj) x(t)  (29)

in which x(t) is the signal observed on all the sensors of the array.The relationship between the variable w=cos(θ−α)cos(Δ) and thecoordinates of the wave vector (u,v) of the equation (2) is as follows:

w=u cos(α)+ν sin(α)  (30)

Knowing w, the incidence Θ becomes a 1D function dependent on theparameter u such that:

$\begin{matrix}{{\underset{\_}{\Theta}(u)} = {\left( {\theta,\Delta} \right) = {{f\left( {u,v} \right)} = {f\left( {u,\frac{w - {u\; {\cos (\alpha)}}}{\sin (\alpha)}} \right)}}}} & (31)\end{matrix}$

in which the function f(u,v) is given by the expression (2). When α=0,the vector of parameter Θ cannot depend on the variable u: In this case,the incidence Θ depends on the variable v with the function Θ(v)=f((w−νsin(α))/ cos(α),ν). In the interests of simplicity in the description ofthe method and without compromising generality, it will be assumed thatit is still possible to write Θ as a function of u. Consequently,

$\begin{matrix}{{\underset{\_}{\Theta}}_{p} = {{{f_{p}\left( u_{p} \right)}\mspace{14mu} {with}\mspace{14mu} {f_{p}(u)}} = {f\left( {u,\frac{w_{p} - {u\; {\cos (\alpha)}}}{\sin (\alpha)}} \right)}}} & (32)\end{matrix}$

With P paths of which at least one group of K_(max) are coherent, theexample described of the method according to the invention contains atleast the following steps:

Step A: Application of a spatial smoothing and/or forward-backwardtechnique to the observation vector x^(l)(t) of the linear array. Aftera 1D goniometry according to the variable w, the incidence parametersw_(p)=cos(θ_(p)−α)cos(Δ_(p)) are obtained for (1≦p≦P). The 1D MUSICcriterion has the following expression:

$\begin{matrix}{{J_{MUSIC}(w)} = \frac{{a^{l}(w)}^{H}\Pi_{b}^{l}{a^{l}(w)}}{{a^{l}(w)}^{H}{a^{l}(w)}}} & (33)\end{matrix}$

in which Π_(b) ^(l) is the noise projector extracted from the smoothedcovariance matrix.

Step B: With K_(max)≦P coherent paths, application of the coherent MUSICmethod described above with the variable Θ={Θ ₁ . . . Θ _(K) _(max) }which is the following function of the variable u={u₁ . . . u_(K) _(max)}

Θ={f ₁(u ₁) . . . f _(K) _(max) (u _(K) _(max) )}  (34)

in which the coherent MUSIC criterion J_(coherent MUSIC) is a functionof the variable u having K_(max) dimensions.

Step C: from K (K being the rank of the covariance matrix R_(x))solutions u _(k) minimizing the function J_(coherent MUSIC) (u) it ispossible to extract the P pairs of incidences (w_(p),u_(p)) for (1≦p≦P)and deduce the incidences (θ_(p),Δ_(p)) therefrom by performing

$\begin{matrix}\left\{ {{\begin{matrix}{\theta_{p} = {{angle}\left( {u_{p} + {j\; v_{p}}} \right)}} \\{\Delta_{p} = {\cos^{- 1}\left( {{u_{p} + {jv}_{p}}} \right)}}\end{matrix}\mspace{14mu} {in}\mspace{14mu} {which}\mspace{14mu} v_{p}} = \frac{w_{p} - {u_{p}{\cos (\alpha)}}}{\sin (\alpha)}} \right. & (35)\end{matrix}$

The preceding steps show that the calculation of a criterion with2K_(max) dimensions for the coherent MUSIC algorithm alone in 2D hasbeen replaced by the calculation of a MUSIC criterion with one dimensionaccording to the parameter w and the cost of calculation of the 1Dcoherent MUSIC criterion with the variable u having _(Kmax) dimensions.The gain in computation power is then equal to

$\begin{matrix}{{Gain} = {\frac{{nb}^{({{2\; K_{\max}} - 1})}}{1 + {nb}^{({K_{\max} - 1})}} \approx {nb}^{K_{\max}}}} & (36)\end{matrix}$

in which nb is the number of points of the meshes of the criteria (MUSICor coherent MUSIC) according to the variables u and v of the componentsof the wave vector. In the general case nb is large while beingproportional to the size of the array (nb>50).

It will be assumed that u _(k) is a solution parameter vector forcoherent MUSIC when

J_(coherent MUSIC)(u _(k))<η(K_(max))  (37)

in which η(K_(max)) is a threshold between 0 and 1 because the criterionJ_(coherent MUSIC) (u) is normalized. When the number K′ of solutions u_(k) is less than the rank K of the covariance matrix R_(x), it can bededuced therefrom that the number of coherent paths is greater thanK_(max). For the case where K′<K_(max) the coherent MUSIC algorithm willbe applied with K_(max)=K_(max)+1. Consequently, the method makes itpossible to jointly estimate the incidences of the paths with the numberof coherent paths.

Similarly, it will be assumed that w_(p) is a solution parameter of the1D goniometry of step A when

J_(MUSIC)(w_(p))<η  (38)

In which η is a threshold between 0 and 1 because the criterionJ_(MUSIC)(w) is normalized.

The method envisages treating the case in which at least two coherentpaths satisfy w_(i)=w_(j) with u_(i)≠u_(j). This problem can be detectedwhen:

the rank of the smoothed covariance matrix remains equal to that of thecovariance matrix Rx;

the MUSIC method does not work on the non-smoothed covariance matrix Rx.

By assuming that there are K_(max) coherent paths and that 1D MUSIC wgives P′<K_(max) coherent path solutions, the method consists incomplementing the incomplete list of P′ elements {w_(p)} with K_(max),−P′ estimation of the initial list of {w_(p)}. More specifically, withK_(max)=2 coherent paths and P′=1 parameter w₁ detected, it is essentialto apply the coherent MUSIC method of step B with w₁=w₁ and w₂=w₁ or theset of parameters {w₁, w₁}. In the case where K_(max)=3 coherent pathsand P′=2 , there are two configurations to which the 1D coherent MUSICstep B must be applied: {w₁, w₂, w₂} and {w₁, w₁, w₂}. Consequently,when P′<K_(max) the step B of the method can be applied several times.There are thus L sets of following incidences w_(p) to which the step Bof the method must be applied:

$\begin{matrix}{{\Omega_{i} = {\Psi\bigcup{X_{i}\mspace{14mu} {for}\mspace{14mu} \left( {1 \leq i \leq L} \right)}}}{{with}\mspace{14mu} \left\{ \begin{matrix}{\Psi = \left\{ {{w_{p}\mspace{14mu} {for}\mspace{14mu} 1} \leq p \leq P^{\prime}} \right\}} \\{{X_{i} \Subset {\Psi \mspace{14mu} {and}\mspace{14mu} {{cardinal}\left( X_{i} \right)}}} = {K_{\max} - P^{\prime}}}\end{matrix} \right.}} & (39)\end{matrix}$

The number L and the sets Ω_(i) can be determined by a conventionalarithmetical process.

The following steps of the method make it possible to estimate thedirection of arrival of P paths in azimuth-elevation bearing in mindthat there is at least one group of K_(max) coherent paths and that thearray is disturbed by mutual coupling of known matrix Z.

Step No. 1: Breakdown the covariance matrix R_(x)=E[x(t) x(t)^(H)] intospecific elements such that

R_(x)=EΛE^(H)

in which E is the matrix with the specific vectors and Λ is a diagonalmatrix consisting of the specific values.

Step No. 2: From the specific values of the matrix Λ, determination ofthe number K of dominant specific values giving the rank of R_(x).

Step No. 3: Following breakdown of the matrix R_(x)

$R_{x} = {{E_{s}\Lambda_{s}E_{s}^{H}} + {E_{b}\Lambda_{b}E_{b}^{H}\mspace{14mu} {with}\mspace{14mu} \left\{ \begin{matrix}{E = \begin{bmatrix}E_{s} & E_{b}\end{bmatrix}} \\{\Lambda = \begin{bmatrix}\Lambda_{s} & 0 \\0 & \Lambda_{b}\end{bmatrix}}\end{matrix} \right.}}$

In which E_(s) and E_(b) are the matrices of the specific vectorsrespectively associated with the signal space bearing in mind that dim(E_(s))=N×K and in which Λ_(s) and Λ_(b) are diagonal matricesrespectively consisting of the specific values of the signal space andthe specific values of the noise space.

Step No 4: Extraction of the non-noise-affected and coupling-freecovariance matrix by performing

$R_{y} = {{Z^{- 1}\left( {R_{x} - {\frac{{trace}\left( \Lambda_{b} \right)}{N - K}I_{N}}} \right)}Z^{{- 1}\; H}}$

In which K is the dimension of the signal space such that K≦P.

Step No. 5: Calculation of the noise projector of the matrix R_(y) inthe following way:

Π_(b) =I _(N) −Z ⁻¹(E _(s)(E _(s) ^(H) Z ^(−1H) Z ⁻¹ E _(s))⁻¹ E _(s)^(H))Z ^(−1H)

Step No. 6: Application of the 2D MUSIC algorithm with the criterionJ_(MUSIC) (Θ)=(

(Θ)^(H Π) _(b) ^(H)

(Θ))/(

(Θ)^(H)

(Θ)) with the vector

(Θ) of the equation (22). Estimation of P₀≦K incidences Θ _(p) (1≦p≦P₀)satisfying J_(MUSIC)(Θ _(p))<η(1). Formation of the set Θ={Θ=Θ ₁ . . . Θ_(P0)} of the non-coherent paths. If P₀<K, go to step No. 7.

Step No. 7: Calculation of the covariance matrix of the linear array byperforming R_(x) ^(l)=P_(roj)R_(x)P_(roj) ^(H) bearing in mind thatx^(l)(t)=P_(roj)x(t).

Step No. 8: Application of one (or both) of the smoothing techniques tothe matrix R_(x) ^(l) of the linear array by performing either {tildeover (R)}_(x)=Σ_(i=1) ^(I)P^(i)R_(x) ^(l)(P^(i))^(H) for smoothing or{tilde over (R)}_(x)=R_(x) ^(l)+Π(R_(x) ^(l))*Π^(T) for theforward-backward.

Step No. 9: From a breakdown into specific elements of the matrix {tildeover (R)}_(x), estimation of the rank P of the signal space and of thenoise projector Π_(b) ^(l)=E_(b)E_(b) ^(H) (step of the MUSIC algorithmreviewed in steps 1 to 3 of this method for the matrix R_(x)).

Step No. 10: Application of the 1D MUSIC algorithm with the criterionJ_(MUSIC) (w)=(a^(l)(w)^(H)Π_(b) ^(l)a^(l)(w))/(a^(l)(w)^(h) a^(l)(w))with the vector a^(l)(w) of the equation (28). Estimation of Pincidences w_(p) (1≦p≦P) satisfying J_(MUSIC)(w_(p))<η.

Step No. 11: Formation of the set Ψ of the incidences w_(p) associatedwith a coherent path such that Ψ={w_(p)≠cos(θ−α)cos(Δ_(i)) in which Θ_(i)={θ_(i), Δ_(i)}∈Θ}

Step No. 12: If P>K, then K_(max)=cardinal(Ψ) and L=1 with Ω₁=Ψ. Go tostep No. 14.

Step No. 13: If P≦K, then K_(max)=K+1 and formation of the L sets ofparameters Ω_(i) of the equation (39) with P′=cardinal(Ψ).

Step No. 14: i=1

Step No. 15: Application of the 1D coherent MUSIC steps B and Cdescribed on page 18 with Θ=f(u)={f₁(u₁) . . . f_(K) _(max) (u_(K)_(max) )}, bearing in mind that f_(p)(u)=f(u,(w_(p)−u cos(α))/ sin(α))in which w_(p)∈Ω_(i). Obtaining of K^(i) incidences Θ _(k) for(1≦k≦K^(i)).

Step No. 16: For k ranging from 1 to K^(i) if Θ _(k)∉Θ then Θ=Θ∪{Θ _(k)}

Step No. 17: i=i+1. If i≦L then return to step No. 14.

One advantage of the method according to the invention is that theminimum number of sensors for estimating the direction of arrival of Kcoherent paths in 2D is lower than with the methods of the prior art,which require a number of sensors greater than 2(K+1), the methodaccording to the invention requiring only a number of sensors greaterthan K.

Another advantage of the method according to the invention is that itmakes it possible to estimate directions of arrival of the paths in 2Dwith larger arrays, which enhances the accuracy of the estimation.

BIBLIOGRAPHY

[1] R O. SCHMIDT, Multiple emitter location and signal parameterestimation, in Proc of the RADC Spectrum Estimation Workshop, GriffithsAir Force Base, New York, 1979, pp. 243-258.

[2] P. Larzabal Application du Maximum de vraisemblance au traitementd'antenne: radio-goniométrie et poursuite de cibles. PhD Thesis,Université de Paris-sud, Orsay, FR, June 1992

[3] B. Ottersten, M. Viberg, P. Stoica and A. Nehorai, Exact and largesample maximum likelihood techniques for parameter estimation anddetection in array processing. In S. Haykin, J. Litva and T J. Shepherseditors, Radar Array Processing, chapter 4, pages 99-151.Springer-Verlag, Berlin 1993.

[4] A. FERREOL, E. BOYER, and P. LARZABAL, <<Low cost algorithm for somebearing estimation in presence of separable nuisance parameters>>,Electronic-Letters, IEE, vol 40, No. 15, pp 966-967, July 2004

[5] S. U. Pillai and B. H. Kwon, Forward/backward spatial smoothingtechniques for coherent signal identification, IEEE Trans. Acoust.,Speech and Signal Processing, vol. 37, pp. 8-15, Jan. 1988

[6] B. Friedlander and A. J. Weiss. Direction Finding using spatialsmoothing with interpolated arrays. IEEE Transactions on Aerospace andElectronic Systems, Vol. 28, No. 2, pp. 574-587, 1992.

[7] A. Ferréol, J.Brugier and Ph.Morgand Method for estimating theangles of arrival of coherent sources by a spatial smoothing techniqueon any array of sensors. Patent published under the number FR 2917180.

1. A method for jointly determining the azimuth angle θ and theelevation angle Δ of the wave vectors of P waves in a system comprisingan array of sensors, a number of waves out of the P waves beingpropagated along coherent or substantially coherent paths between asource and said sensors, the method comprising at least the followingsteps: selecting a subset of sensors from said sensors to form a linearsubarray (1002) of sensors; applying, to the signals from the chosensubarray, an algorithm according to a single dimension to decorrelatethe sources of the P waves; determining a first component w of said wavevectors by applying, to the signals observed on the sensors of thechosen subarray, a goniometry algorithm according to the singledimension w; determining a second component u of said wave vectors byapplying a goniometry algorithm according to the single dimension u tothe signals from the entire array of sensors; determining θ and Δ from wand u.
 2. The method as claimed in claim 1, wherein: the sensors formingthe subarray are chosen such that at least a portion of the subarray isunchanging by translation; a spatial smoothing algorithm is applied todecorrelate the sources of the P waves.
 3. The method as claimed inclaim 1, wherein: the sensors forming the subarray are chosen such thatat least a portion of the subarray includes a center of symmetry; aforward-backward algorithm is applied to decorrelate the sources of theP waves.
 4. The method as claimed in claim 1 wherein the determinedfirst component w of the wave vectors is the projection, on the axisformed by the linear subarray, of the projection of the wave vectors onthe plane formed by the array of sensors.
 5. The method as claimed inclaim 1 further comprising at least the following steps: calculating thecovariance matrix R_(x) on the entire array of sensors; extracting fromR_(x), the covariance matrix R_(x)′ corresponding to the chosen linearsubarray; applying a source decorrelation algorithm to R_(x)′;estimating, for each path p, the values of the first componentw_(p)=cos(θ_(p)−α).cos(Δ_(p)) by applying a 1D goniometry algorithm tothe decorrelated matrix R_(x)′, α being the azimuth orientation angle ofthe axis formed by the linear subarray; estimating the values of thesecond component u_(p)=cos(θ_(p)).cos(Δ_(p)), for each path p, byapplying a 1D goniometry algorithm to the matrix R_(x); determining,from the values of the pairs (w_(p), u_(p)), the values of theazimuth-elevation pairs (θ_(p), Δ_(p)).
 6. The method as claimed inclaim 1 wherein the goniometry algorithm used to determine the firstcomponent w of each wave vector P is the MUSIC algorithm, the criterionJ_(MUSIC) to be minimized to determine said component w being equal to${{J_{MUSIC}(w)} = \frac{{a^{l}(w)}^{H}\Pi_{b}^{l}{a^{l}(w)}}{{a^{l}(w)}^{H}{a^{l}(w)}}},$in which Π_(b) ^(l) is the noise projector extracted from thedecorrelated covariance matrix R_(x)′ corresponding to the chosen linearsubarray, and a(w)^(l) represents the response of the chosen subarray tothe incident waves P.
 7. The method as claimed in claim 1 wherein thegoniometry algorithm used to determine the second component is thecoherent MUSIC algorithm in a single dimension, the criterion to beminimized being:${J_{{MUSICcoherent}\mspace{14mu} {MUSIC}}\left( \underset{\_}{\Theta} \right)} = \frac{\det \left( {{D\left( \underset{\_}{\Theta} \right)}^{H}\Pi_{b}{D\left( \underset{\_}{\Theta} \right)}} \right)}{\det \left( {{D\left( \underset{\_}{\Theta} \right)}^{H}{D\left( \underset{\_}{\Theta} \right)}} \right)}$in which Θ={f_(l)(u_(l)) . . . f_(K) _(max) (u_(K) _(max) )}, with${{f_{p}(u)} = {f\left( {u,\frac{w_{p} - {u\; {\cos (\alpha)}}}{\sin (\alpha)}} \right)}},$α being the azimuth orientation angle of the axis formed by the linearsubarray, D(Θ) being a vector equal to [a(Θ ₁) . . . a(Θ _(Kmax)], a(Θ₁) being the response of the array of sensors to the path of index i,K_(max) being the maximum number of coherent paths.
 8. The method asclaimed in claim 1 wherein the goniometry algorithm used to determinethe second component is the maximum likelihood algorithm.
 9. The methodas claimed in claim 5, the array being disturbed by mutual coupling ofknown matrix Z, further comprising a step for eliminating the couplingexecuted prior to the steps for estimating the values of the componentsw and u, said step for eliminating the coupling determining a covariancematrix that is cleaned of noise by applying the following processing tothe covariance matrix: Z⁻¹(R_(x)−σ²|)Z^(−1H), σ² being the estimatednoise level.
 10. The method as claimed in claim 5, wherein thedetermination of the pairs of values (θ_(p), Δ_(p)) from the values ofpairs (w_(p), u_(p)) is performed as follows:$\quad\left\{ \begin{matrix}{\theta_{p} = {{angle}\left( {u_{p} + {j\; v_{p}}} \right)}} \\{\Delta_{p} = {\cos^{- 1}\left( {{u_{p} + {jv}_{p}}} \right)}}\end{matrix}\mspace{11mu} \right.$ in which$v_{p} = {\frac{w_{p} - {u_{p}{\cos (\alpha)}}}{\sin (\alpha)}.}$